OHI British Columbia | OHI Science | Citation policy
knitr::opts_chunk$set(fig.width = 6, fig.height = 4, fig.path = 'figs/',
echo = TRUE, message = FALSE, warning = FALSE)
library(raster)
library(data.table)
library(sf)
source('https://raw.githubusercontent.com/oharac/src/master/R/common.R') ###
### includes library(tidyverse); library(stringr); dir_M points to ohi directory
dir_git <- '~/github/spp_risk_dists'
### goal specific folders and info
dir_data <- file.path(dir_git, 'data')
dir_o_anx <- file.path(dir_O, 'git-annex/spp_risk_dists')
source(file.path(dir_git, 'setup/common_fxns.R'))
### Gall-Peters doesn't have an EPSG?
gp_proj4 <- '+proj=cea +lon_0=0 +lat_ts=45 +x_0=0 +y_0=0 +ellps=WGS84 +units=m +no_defs'Create a set of maps of the distribution of biodiversity intactness - all species assessed and mapped by IUCN. These maps are generated at 10 km2 resolution in a Gall-Peters projection. These maps will be generated using all available species:
A selection of these maps will be generated for taxonomic groups and range sizes in a separate Rmd.
Future iterations may include:
IUCN Red List spatial data download IUCN Red List API Gina Ralph (IUCN)
From the 1a and 1b biodiversity risk map scripts, gather the rasters for the various maps:
risk_un_rast <- raster(file.path(dir_git, 'output', 'mean_risk_raster_comp.tif'))
risk_rr_rast <- raster(file.path(dir_git, 'output', 'mean_rr_risk_raster_comp.tif'))
threat_un_rast <- raster(file.path(dir_git, 'output', 'pct_threat_raster_comp.tif'))
threat_rr_rast <- raster(file.path(dir_git, 'output', 'sr_rr_pct_threat_raster_comp.tif'))
n_spp_rast <- raster(file.path(dir_git, 'output', 'n_spp_risk_raster_comp.tif'))
n_rr_rast <- raster(file.path(dir_git, 'output', 'sr_rr_risk_raster_comp.tif'))four_panel <- function(map1_rast, map2_rast,
limits = c(0, 1),
colors, values,
labels, breaks,
plot_labs = c('A', 'B')) {
land_poly <- sf::read_sf(file.path(dir_git, 'spatial/ne_10m_land/ne_10m_land.shp')) %>%
st_transform(gp_proj4)
map_df <- data.frame(val_1 = values(map1_rast),
val_2 = values(map2_rast)) %>%
cbind(coordinates(map1_rast)) %>%
filter(!is.na(val_1))
map1 <- ggplot(map_df) +
geom_raster(aes(x, y, fill = val_1), show.legend = FALSE) +
geom_sf(data = land_poly, aes(geometry = geometry),
fill = 'grey96', color = 'grey40', size = .10) +
ggtheme_map() +
theme(plot.margin = unit(c(.2, 0, .1, .5), units = 'cm')) +
coord_sf(datum = NA) + ### ditch graticules
scale_fill_gradientn(colors = colors, values = values, limits = limits,
labels = labels, breaks = breaks,
guide = guide_colourbar(label.position = 'left',
label.hjust = 1)) +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0))
map2 <- ggplot(map_df) +
geom_raster(aes(x, y, fill = val_2), show.legend = FALSE) +
geom_sf(data = land_poly, aes(geometry = geometry),
fill = 'grey96', color = 'grey40', size = .10) +
ggtheme_map() +
theme(plot.margin = unit(c(.1, 0, .2, .5), units = 'cm')) +
coord_sf(datum = NA) + ### ditch graticules
scale_fill_gradientn(colors = colors, values = values, limits = limits,
labels = labels, breaks = breaks,
guide = guide_colourbar(label.position = 'left',
label.hjust = 1)) +
scale_x_continuous(expand = c(0, 0)) +
scale_y_continuous(expand = c(0, 0))
###################
### set up a dataframe of values to craft a color bar using geom_segment
colorbar_df <- data.frame(x = seq(0, 1, .001), y = -1)
dens1 <- ggplot(map_df) +
ggtheme_plot() +
geom_vline(xintercept = mean(map_df$val_1, na.rm = TRUE)) +
geom_segment(data = colorbar_df,
aes(x = x, xend = x, color = x),
y = 0, yend = -10, size = 1,
show.legend = FALSE) +
geom_density(aes(x = val_1, ..scaled..), alpha = .3, size = .25, fill = 'grey30') +
scale_color_gradientn(colors = colors, values = values, limits = limits,
labels = labels, breaks = breaks) +
theme(axis.text.x = element_blank(),
axis.title = element_blank(),
panel.grid.major.x = element_blank(),
plot.margin = unit(c(1, 0, 1, 0), units = 'cm')) +
scale_x_continuous(labels = labels, breaks = breaks, limits = limits,
expand = c(0, 0)) +
scale_y_continuous(limits = c(-.4, 1)) +
coord_flip()
dens2 <- ggplot(map_df) +
ggtheme_plot() +
geom_vline(xintercept = mean(map_df$val_2, na.rm = TRUE)) +
geom_segment(data = colorbar_df,
aes(x = x, xend = x, color = x),
y = 0, yend = -10, size = 1,
show.legend = FALSE) +
geom_density(aes(x = val_2, ..scaled..), alpha = .3, size = .25, fill = 'grey30') +
scale_color_gradientn(colors = colors, values = values, limits = limits,
labels = labels, breaks = breaks) +
theme(axis.text.x = element_blank(),
axis.title = element_blank(),
panel.grid.major.x = element_blank(),
plot.margin = unit(c(1, 0, 1, 0), units = 'cm')) +
scale_x_continuous(labels = labels, breaks = breaks, limits = limits,
expand = c(0, 0)) +
scale_y_continuous(limits = c(-.4, 1)) +
coord_flip()
panel_top <- cowplot::plot_grid(map1, dens1,
axis = 'b',
rel_widths = c(5, 1))
panel_btm <- cowplot::plot_grid(map2, dens2,
axis = 'b',
rel_widths = c(5, 1))
four_panel <- cowplot::plot_grid(panel_top, panel_btm,
labels = plot_labs,
nrow = 2, ncol = 1, align = 'v')
}### aggregate rasters for faster testing:
# risk_un_rast <- risk_un_rast %>% aggregate(10)
# risk_rr_rast <- risk_rr_rast %>% aggregate(10)
mean_four_panel <- four_panel(risk_un_rast, risk_rr_rast,
colors = risk_cols,
values = risk_vals,
labels = risk_lbls,
breaks = risk_brks,
plot_labs = c('A', 'B'))### aggregate rasters for faster testing:
# threat_un_rast <- threat_un_rast %>% aggregate(10)
# threat_rr_rast <- threat_rr_rast %>% aggregate(10)
threat_four_panel <- four_panel(threat_un_rast, threat_rr_rast,
colors = thr_cols,
values = thr_vals,
breaks = thr_brks,
labels = thr_lbls,
plot_labs = c('C', 'D'))four_maps <- cowplot::plot_grid(mean_four_panel, threat_four_panel,
nrow = 1, ncol = 2, align = 'v')
fig1_path <- file.path(dir_git, 'ms_figures/fig1_comp_assessed.png')
ggsave(plot = four_maps,
filename = fig1_path,
width = 10, height = 5, units = 'in', dpi = 300)From the 1a and 1b biodiversity risk map scripts, gather the rasters for the various maps:
risk_un_rast <- raster(file.path(dir_git, 'output',
'mean_risk_raster_all.tif'))
risk_rr_rast <- raster(file.path(dir_git, 'output',
'mean_rr_risk_raster_all.tif'))
threat_un_rast <- raster(file.path(dir_git, 'output',
'pct_threat_raster_all.tif'))
threat_rr_rast <- raster(file.path(dir_git, 'output',
'sr_rr_pct_threat_raster_all.tif'))
n_spp_rast <- raster(file.path(dir_git, 'output',
'n_spp_risk_raster_all.tif'))
n_rr_rast <- raster(file.path(dir_git, 'output',
'sr_rr_risk_raster_all.tif'))### aggregate rasters for faster testing:
# risk_un_rast <- risk_un_rast %>% aggregate(10)
# risk_rr_rast <- risk_rr_rast %>% aggregate(10)
mean_four_panel <- four_panel(risk_un_rast, risk_rr_rast,
colors = risk_cols,
values = risk_vals,
labels = risk_lbls,
breaks = risk_brks,
plot_labs = c('A', 'B'))### aggregate rasters for faster testing:
# threat_un_rast <- threat_un_rast %>% aggregate(10)
# threat_rr_rast <- threat_rr_rast %>% aggregate(10)
threat_four_panel <- four_panel(threat_un_rast, threat_rr_rast,
colors = thr_cols,
values = thr_vals,
breaks = thr_brks,
labels = thr_lbls,
plot_labs = c('C', 'D'))four_maps <- cowplot::plot_grid(mean_four_panel, threat_four_panel,
nrow = 1, ncol = 2, align = 'v')
fig1all_path <- file.path(dir_git, 'ms_figures/figSI_all_species.png')
ggsave(plot = four_maps,
filename = fig1all_path,
width = 10, height = 5, units = 'in', dpi = 300)